package com.power.sifa.util;

import java.util.List;

/**
 * @ClassName GeoUtil
 * @Author xiongzw
 * @Date 2023/12/10 13:48
 * @Description
 * @Version: 1.0
 **/

public class GeoUtil {

    /**
     * 地球半径（单位：m）
     */
    private static final Double Zc = 6378136.49;

    public static Double toRadians(Double angdeg) {
        return angdeg * (Math.PI / 180);
    }

    public static Double Yd (Double angdeg) {
        return angdeg / (Math.PI / 180);
    }

    public static Double Hf (Double a, Integer b, Integer c) {
        return a >= b && a <= c ? a : ((a - b) % (c - b) + (c - b)) % (c - b) + b;
    }

    /**
     * 两个坐标点间的距离
     * @param lnga 起始点经度
     * @param lata 起始点纬度
     * @param lngb 终点经度
     * @param latb 终点纬度
     * @return distance （单位：米）
     */
    public static Double computeDistanceBetween(Double lnga ,Double lata ,Double lngb ,Double latb) {
        Double c = toRadians(lata),d=toRadians(latb);
        return (2*Math.asin(Math.sqrt(Math.pow(Math.sin((c-d)/2),2)+Math.cos(c)*Math.cos(d)*Math.pow(Math.sin((toRadians(lnga)-toRadians(lngb))/2),2))))*Zc;
    }

    /**
     * 从一个坐标到另一个坐标的航向。航向是指从一个坐标指向另一个坐标的向量与正北方向的夹角，范围为[-180,180)。
     * @param lnga 起始点经度
     * @param lata 起始点纬度
     * @param lngb 终点经度
     * @param latb 终点纬度
     * @return angle
     */
    public static Double computeHeading(Double lnga ,Double lata ,Double lngb ,Double latb) {
        Double c = toRadians(lata),d=toRadians(latb),e=toRadians(lngb)-toRadians(lnga);
        return Hf(Yd(Math.atan2(Math.sin(e)*Math.cos(d),Math.cos(c)*Math.sin(d)-Math.sin(c)*Math.cos(d)*Math.cos(e))),-180,180);
    }

    /**
     * 通过起始点坐标、距离以及航向算出终点坐标
     * @param lng 起始点经度
     * @param lat 起始点纬度
     * @param distance 距离（单位：m）
     * @param heading 航向角度 [-180,180]
     * @return [lng, lat]
     */
    public static Double[] computeOffset (Double lng ,Double lat, Double distance, Double heading) {
        distance /= Zc;
        heading = toRadians(heading);
        Double e = toRadians(lat);
        Double d = Math.cos(distance);
        distance = Math.sin(distance);
        Double f = Math.sin(e);
        e=Math.cos(e);
        Double g = d*f+distance*e*Math.cos(heading);
        return new Double[]{Yd(toRadians(lng)+Math.atan2(distance*e*Math.sin(heading),d-f*g)), Yd(Math.asin(g))};
    }

    /**
     * 通过终点坐标、距离以及航向算出起始点坐标
     * @param lng 终点经度
     * @param lat 终点纬度
     * @param distance 距离（单位：m）
     * @param heading 航向角度 [-180,180]
     * @return
     */
    public static Double[] computeOffsetOrigin (Double lng ,Double lat, Double distance, Double heading) {
        heading = toRadians(heading);
        distance /= Zc;
        Double d=Math.cos(distance);
        Double e = Math.sin(distance)*Math.cos(heading);
        distance = Math.sin(distance)*Math.sin(heading);
        heading=Math.sin(toRadians(lat));
        Double f = e*e*d*d+d*d*d*d-d*d*heading*heading;
        if(0>f) {
            return null;
        }
        Double g= e*heading+Math.sqrt(f);
        g=g/(d*d+e*e);
        Double h=(heading-e*g)/d;
        g=Math.atan2(h,g);
        if(g<-Math.PI/2||g>Math.PI/2) {
            g=e*heading-Math.sqrt(f);
            g=Math.atan2(h,g/(d*d+e*e));
        }
        return g<-Math.PI/2||g>Math.PI/2?null:new Double[]{Yd(toRadians(lng)-Math.atan2(distance,d*Math.cos(g)-e*Math.sin(g))), Yd(g)};
    }

    /**
     * 获取坐标中心点
     * @param pointList List<Double[lng,lat]>
     * @return
     */
    public static Double[] computeCenter(List<Double[]> pointList){
        Double[] center= null;
        if (pointList == null || pointList.isEmpty()){
            return center;
        }
        center = new Double[2];

        Double lng = 0.0;
        Double lat = 0.0;
        for(Double[] point : pointList){
            if (point[0] != null && point[1] != null) {
                lng += (Double)point[0];
                lat += (Double)point[1];
            }
        }
        center[0]= lng/pointList.size();
        center[1]= lat/pointList.size();
        return center;
    }

    /**
     * 计算两点间的球面距离
     * @param alng
     * @param alat
     * @param blng
     * @param blat
     * @return
     */
    public static Double calcP2P(double alng, double alat, double blng, double blat) {
        // 判断点是在北半球还是南半球，本程序中输入的数据若为负则表示在南边球
        double distance = 0.0D;

        double _alat = toRadians(alat);
        double _alng = toRadians(alng);
        double _blat = toRadians(blat);
        double _blng = toRadians(blng);

        double c = Math.sin(_alat) * Math.sin(_blat) + Math.cos(_alat) * Math.cos(_blat) * Math.cos(_alng - _blng); // Java中三角函数角度以弧度制表示
        if (c > 1) {
            c = 1;
        }
        distance = Math.acos(c) * Zc; // 弧长公式：弧长 = 弧度 * 半径
        if (distance <= 0.01) { // GPS误差
            distance = 0.0D;
        }
        return distance;
    }

    /**
     * 计算点到直线的最短距离（单位：米）
     * @param x 点经度
     * @param y 点纬度
     * @param x1 线起点经度
     * @param y1 线起点纬度
     * @param x2 线终点经度
     * @param y2 线终点纬度
     * @return
     */
    public static Double computeDistanceToLine(double x, double y, double x1, double y1, double x2, double y2) {
        double d1 = (x2 - x1) * (x - x1) + (y2 - y1) * (y - y1);
        if (d1 < 0) {
            return calcP2P(x, y, x1, y1);
        }
        double d2 = (x2 - x1) * (x2 - x1) + (y2 - y1) * (y2 - y1);
        if (d1 >= d2) {
            return calcP2P(x, y, x2, y2);
        }
        double r = d1 / d2;
        double px = x1 + (x2 - x1) * r;
        double py = y1 + (y2 - y1) * r;
        return calcP2P(x, y, px, py);
    }

    /**
     * 判断点是否在多边形内（只支持简单多边形， 即一维数组）
     * @param lng 点经度
     * @param lat 点纬度
     * @param polygon 多边形
     * @return
     */
    public static boolean isPointInPolygon(double lng, double lat, List<double[]> polygon) {
        if (polygon == null || polygon.isEmpty() || polygon.size()<3){
            return false;
        }
        List<double[]> pts = polygon;
        int N = pts.size();
        boolean boundOrVertex = true;
        int intersectCount = 0;
        double precision = 2e-10;
        double[] p1, p2;
        double[] p = new double[]{lng, lat};
        p1 = pts.get(0);
        for (int i = 1; i <= N; ++i) {
            if (p[1] == p1[1] || p[0] == p1[0]) {
                return boundOrVertex;
            }
            p2 = pts.get(i % N);
            if (p[1] < Math.min(p1[1], p2[1]) || p[1] > Math.max(p1[1], p2[1])) {
                p1 = p2;
                continue;
            }
            if (p[1] > Math.min(p1[1], p2[1]) && p[1] < Math.max(p1[1], p2[1])) {
                if (p[0] <= Math.max(p1[0], p2[0])) {
                    if (p1[1] == p2[1] && p[0] >= Math.min(p1[0], p2[0])) {
                        return boundOrVertex;
                    }
                    if (p1[0] == p2[0]) {
                        if (p1[0] == p[0]) {
                            return boundOrVertex;
                        } else {
                            ++intersectCount;
                        }
                    } else {
                        double xinters = (p[1] - p1[1]) * (p2[0] - p1[0]) / (p2[1] - p1[1]) + p1[0];
                        if (Math.abs(p[0] - xinters) < precision) {
                            return boundOrVertex;
                        }
                        if (p[0] < xinters) {
                            ++intersectCount;
                        }
                    }
                }
            } else {
                if (p[1] == p2[1] && p[0] <= p2[0]) {
                    double[] p3 = pts.get((i + 1) % N);
                    if (p[1] >= Math.min(p1[1], p3[1]) && p[1] <= Math.max(p1[1], p3[1])) {
                        ++intersectCount;
                    } else {
                        intersectCount += 2;
                    }
                }
            }
            p1 = p2;
        }
        if (intersectCount % 2 == 0) {
            return false;
        } else {
            return true;
        }
    }

    /**
     * 计算多边形面积（平方米）
     * @param points
     * @return
     */
    public static double computeArea(List<double[]> polygon) {
        if (polygon == null || polygon.isEmpty() || 3>polygon.size()){
            return 0.0;
        }
        Double c = Zc * (Math.PI/180);
        Double d=0.0;
        int e = polygon.size();
        int x = 0;
        for(int f=0; f< e-1; f+=1){
            double[] h = polygon.get(f);
            double[] k = polygon.get(f+1);
            double u = h[0] * c * Math.cos(h[1] * (Math.PI / 180));
            double g = h[1] * c;
            double v = k[0] * c * Math.cos(k[1] * (Math.PI / 180));
            d = d + (u * k[1] * c - v * g);
            x = f;
        }
        double[] fx = polygon.get(x+1);
        double[] bx = polygon.get(0);
        Double ex = fx[0] * c * Math.cos(fx[1] * (Math.PI / 180));
        Double f = fx[1] * c;
        Double k = bx[0] * c * Math.cos(bx[1] * (Math.PI / 180));
        d += ex * bx[1] * c - k * f;
        return 0.5 * Math.abs(d);
    }

}
